Segregation by membrane rigidity in flowing binary suspensions of elastic capsules 
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Spatial segregation in the wall normal direction is investigated in suspensions containing a binary 
mixture of Neo-Hookean capsules subjected to pressure driven flow in a planar slit. The two compo- 
nents of the binary mixture have unequal membrane rigidities. The problem is studied numerically 
using an accelerated implementation of the boundary integral method. The effect of a variety of 
parameters was investigated, including the capillary number, rigidity ratio between the two species, 
volume fraction, confinement ratio, and the number fraction of the more floppy particle ^/ in the 
mixture. It was observed that in suspensions of pure species, the mean wall normal positions of the 
stiff and the floppy particles are comparable. In mixtures, however, the stiff particles were found to 
be increasingly displaced toward the walls with increasing Xf, while the floppy particles were found 
to increasingly accumulate near the centerline with decreasing Xf. This segregation behavior was 
universally observed independent of the parameters. The origin of this segregation is traced to the 
effect of the number fraction X/ on the localization of the stiff and the floppy particles in the near 
wall region - the probability of escape of a stiff particle from the near wall region to the interior is 
greatly reduced with increasing Xf, while the exact opposite trend is observed for a floppy particle 
with decreasing Xf. Simple model studies on heterogeneous pair collisions involving a stiff and a 
floppy particle mechanistically explain the contrasting effect of X/ on the near wall localization of 
the two species. The key observation in these studies is that the stiff particle experiences much larger 
cross-stream displacement in heterogeneous collisions than the floppy particle. A unifled mechanism 
incorporating the wall- induced migration of deformable particles away from the wall and the particle 
fluxes associated with heterogeneous and homogeneous pair collisions is presented. 



I. INTRODUCTION 

Flowing suspensions of mixtures of particles with dif- 
ferent rigidity are naturally encountered in blood flow. 
Blood is primarily a suspension of red blood cells (RBCs) 
with small amounts of white blood cells (WBCs) and 
platelets. Normal RBCs are highly deformable, and can 
thus pass easily through small capillaries, while leuko- 
cytes (WBCs) and platelets are typically much stiffer. 
The leukocytes and platelets are typically found in en- 
hanced concentration near the walls - this property is 
commonly known as margination and is critical for the 
physiological responses of inflammation and hemostasis 

[2] . It is believed that the rigidity of the leukocytes 
and platelets play an important role in this margination 
phenomena. In contrast to the benefits arising from the 
stiff nature of the leukocytes and platelets, the stiffening 
of RBCs in various diseased states like malaria ^ and 
sickle cell disease 3 leads to various deleterious effects. 
A particular lethal complication of malaria is cerebral 
malaria, in which infected RBCs are found to sequester 
in the brain microcirculation Q. This sequestration re- 
sults from the cytoadherence of the infected RBCs to the 
vessel walls of the post capillary venules [5] , and this ef- 
fect is only expected to be exacerbated by the altered 
margination properties of the stiffened RBCs. Recently, 
there has also been a particular emphasis on drug deliv- 
ery with particles via the bloodstream. A common target 
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of such particles is the vascular walls f6j, so in this con- 
text an ideal drug delivery particle must have an inherent 
propensity to segregate toward the walls. As such the de- 
sign of an efficient drug delivery particle (via shape, size, 
and deformability) is an interesting problem in itself. 

The segregation properties of leukocytes, platelets, and 
infected RBCs can also be employed for their separation 
or detection in biomimetic microfluidic devices. For ex- 
ample, Shevkoplyas et al. [3] developed a network of rect- 
angular microfiuidic devices for the separation of leuko- 
cytes from whole blood. The process of selectively re- 
moval of leukocytes is known as leukapheresis and has 
potential applications in decreasing the count of WBCs 
in patients with leukemia or for the removal of WBCs 
for transfusion Q. It was noted above that RBCs in- 
fected with the malaria parasite become stiff. Hou et al. 
employed this property to separate malaria-infected 
RBCs from the rest of the blood via the difference in 
margination properties based on stiffness. They showed 
that even a small difference in deformability is sufficient 
to marginate the stiff infected RBCs, thereby enabling 
their separation. A difference between that study and 
those on leukocytes is that not only are leukocytes stiffer 
than RBCs, they are also larger, which might play a role 
in their margination [lo| . The same is true for margina- 
tion studies on platelets, which are stiffer and smaller 
than RBCs 

There have been several previous simulation-based ef- 
forts to investigate the mechanism(s) of the margina- 
tion behavior of leukocytes and platelets. Munn and 
Dupin [l^ studied leukocyte margination and concluded 
that the formation of stacks of discoidal platelets up- 
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stream of the leukocytes leads to the margination be- 
havior. Freund [lOj saw a similar behavior in his study 
and suggested that the size difference between the RBCs 
and the leukocytes could be playing a role. Munn and 
Dupin [ill also studied suspensions of mixtures of parti- 
cles with different rigidity and found the stiffer particles 
to marginate toward the walls. They explained the phe- 
nomena based on single particle arguments, namely that 
the more deformable particle has a greater tendency to 
migrate toward the centerline and consequently displaces 
the stiffer particles. However, in a suspension, it is ex- 
pected that the particle-particle collisions will play an 
important, perhaps dominant, role, and hence that sin- 
gle particle arguments may not be sufficient to provide 
a full explanation of the phenomenon. Crowl and Fogel- 
son 12] studied numerically the margination behavior of 
platelets in 2-D suspensions of blood. They described 
the phenomenon with a drift-diffusion model, which has 
been used previously for interpreting experimental results 
on the margination behavior of platelets To obtain 
agreement b etween the model an d the detailed simula- 
tion results, ICrowl and FogelsonI had to assume a drift 
velocity of the platelets toward the wall when they are 
close to walls in the cell free layer. They attributed this 
to the one-sided collision of platelets with RBCs, as no 
RBCs are present on the other side (on the side of the 
wall) to balance this drift. There is, of course, nothing 
special about the platelets in the near wall region, as any 
particle in this region will undergo one-sided collisions, 
and as a result get pushed toward the wall. This drift 
term, therefore, does not provide any significant mech- 
anistic insight into the segregation behavior, certainly 
not for segregation in suspensions of equal sized parti- 
cles. In a similar study, using numerical simulations, 
Zhao and Shaqfeh ll| investigated the margination of 
platelets in blood flow. They concluded that the velocity 
fluctuations are responsible for this behavior - once the 
platelets marginate toward the walls, the small velocity 
fluctuations there imply that they can't return back to 
flow. This argument, again, rests on the small size of the 
platelets due to which they are influenced only by the 
local velocity fluctuations, unlike the much larger RBCs. 
Consequently, the segregation between equal sized par- 
ticles based on rigidity, if any, cannot be appropriately 
explained by this mechanism. 

The discussion above clearly establishes the impor- 
tance of rigidity based flow induced segregation in mix- 
tures of particles. At the same time, it is also obvious 
that this phenomena is very poorly understood. In most 
previous works, the key aspects of the proposed mech- 
anisms appear to rest on the shape and the size dis- 
crepancy between the mixture constituents - the role of 
rigidity is not entirely clear and more broadly, no gen- 
eral formalism has been proposed for understanding the 
segregation phenomena. In the present effort, we seek 
to delineate the effect of rigidity on flow induced seg- 
regation by focusing on binary mixtures of deformable 
capsules with unequal rigidities (i.e. stiff and floppy par- 



ticles), but with identical size and spherical rest shape. 
The flow problem involving suspensions of deformable 
capsules is investigated numerically using an accelerated 
implementation of the boundary integral method. The 
specific system geometry and bulk flow considered here 
is pressure driven flow in a planar slit. The detailed sim- 
ulation results indicate two key observations: (i) with in- 
creasing fraction of floppy particles in the suspension, the 
stiff particles get increasingly localized in the near wall 
region, and (ii) with increasing fraction of stiff particles 
in the suspension, the floppy particles are increasingly 
depleted in the near wall region. As a result of this be- 
havior, stiff particles segregate in the near wall region in 
a suspension of primarily floppy particles, while floppy 
particles segregate around the centerline in a suspension 
of primarily stiff particles. The degree of localization of 
stiff and floppy particles in the near wall region as a func- 
tion of the fraction of the other particle type is found to 
be greatly influenced by the cross-stream displacement of 
that species in heterogenous pair collisions, i.e. collisions 
involving a stiff and a floppy particle. In particular, in 
heterogenous pair collisions, the stiff particle undergoes 
much larger cross-stream displacement than the floppy 
particle. As a result, a stiff particle in the near wall re- 
gion gets increasingly localized as the fraction of floppy 
particles in the suspension increases, while the reverse 
is true for a floppy particle in the near wall region as 
the fraction of stiff particles in the suspension increases. 
All these observations can be successfully explained with 
a unified mechanism that incorporates the wall-induced 
migration of deformable suspended particles away from 
the wall and the particle fluxes associated with heteroge- 
neous and homogeneous pair collisions. 

The organization of this article is as follows. In Sec. 
we formulate the problem and discuss the numeri- 
cal solution procedure. Next, in Sec. (jllip . we present 
detailed results for one of the parameter sets explored 
in this work. This is followed by a detailed analysis of 
the segregation mechanism in Sec. pV[) . Concluding re- 
marks are presented in Sec. (|V]). Results for several other 
parameter sets are presented in the appendix. 



II. PROBLEM FORMULATION AND 
IMPLEMENTATION 



Fluid velocity calculation: Boundary integral 
method 



We consider a suspension of fluid-filled deformable 
Neo-Hookean capsules between two parallel plates as 
shown in Fig. ([1]). Both the suspending fiuid and the 
fiuid inside the capsules are assumed to be Newtonian 
and incompressible with the same viscosity, fj,, i.e. the 
viscosity ratio A is unity. We further assume that the 
Reynolds number for the problem is sufficiently small 
that the fluid motion is governed by the Stokes equa- 
tion. Under these assumptions, one may write the fluid 
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following expressions: 
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FIG. 1: (color online) A schematic of the problem set up. 



velocity u at any point xg in the domain as fl5|: 



1 ^ f 

Mj(xo) = (xo) - - — V / A/i(x)Gji(xo,x)dS'(x) 

(1) 

where u°°(xo) is the undisturbed pressure driven velocity 
at point xq, S"" denotes the surface of particle n, Af(x) 
is the hydrodynamic traction jump across the interface, 
while G is the Green's function for the Stokes equation 
in the geometry of interest. Note that the sum in the 
above expression is over all the Np particles in the sys- 
tem. A crucial aspect of the above formulation is that 
the Green's function G is taken to satisfy the boundary 
conditions imposed at the system boundaries, so the inte- 
grals above only involve the internal (interfacial) bound- 
aries; if the Green's function for any other geometry is 
employed (e.g. free-space), additional integrals over the 
domain boundaries will arise in Eq. ([1]). For the present 
case of a slit geometry, we have periodic boundary con- 
ditions in X and z directions with spatial periods and 
Lz, respectively (Fig. [T|). In the y direction, we have 
no slip velocity boundary conditions at the two walls at 
y = Oandy^H (Fig. [J). 

Traditional implementations of the boundary integral 
method typically have a scaling of at least 0{N^), where 
N is proportional to the product of number of particles 
Np and the number of elements A^a employed to dis- 
cretize the particle surface. In the present effort, we em- 
ploy an accelerated implementation of the boundary in- 
tegral method with a computational cost of 0(A^log A^) 
for the slit geometry [IB] • The acceleration in our imple- 
mentation is provided by the use of the General Geom- 
etry Ewald like method (GGEM) [17i], which essentially 
provides fast solution for the velocity and the pressure 
fields driven by a collection of point forces in an arbi- 
trary geometry - finding these velocity and pressure fields 
is at the heart of the solution procedure of the bound- 
ary integral method [l^ . The key aspect of the GGEM 
methodology is to split a Dirac-delta force density into a 
smooth quasi-Gaussian global density Pg(r) and a second 
local density Pi(r); these are respectively given by the 



pi{r) ^d{r)~pg{r), 



(2a) 



(2b) 



where represents a length scale over which the delta- 
function density 6(r) has been smeared using the quasi- 
Gaussian form above, while r is position vector relative 
to the pole of the singularity. It is important to em- 
phasize that the total density remains a J-function, i.e. 
Pg{r) + pi{r) = 6{r). The solution associated with a lo- 
cal density, which is known analytically, is short ranged, 
and is neglected beyond a length scale of 0{a~^) from 
its pole. We note that the local solution is obtained 
assuming free-space boundary conditions. The solution 
associated with the global density is numerically com- 
puted, while ensuring that the boundary conditions as- 
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sociated with the overall problem are satisfied 
Complete details of this method can be found in 
in short, for the slit geometry of interest here we have 
implemented a spectral Galerkin scheme with Fourier se- 
ries representation in the periodic x and z directions and 
Chebyshev polynomial series representation in the y di- 
rection; the number of corresponding modes employed in 
X, y, and z directions are denoted by N^, Ny, and A^z, 
respectively. An important parameter controlling the er- 
ror of the GGEM solution is ahm , where hm, is the mean 
mesh spacing associated with the global solution proce- 
dure (e.g. hm = Lx/Nx). Based on extensive tests in 
(l6j . we have set ahm = 0.5 in this study with the mean 
mesh spacing being equal in all three directions. The 
value of a was held constant in all cases: a — 4/ a, where 
a is the radius of the spherical capsules at rest. 



B. Surface discretization and membrane mechanics 

The surface of each of the particles is discretized into 
A^A triangular elements with linear basis functions em- 
ployed over each element. In the present study, we took 
A^A = 320. As will be seen later, long simulations are re- 
quired to reach a steady state, and this mesh resolution 
keeps the computational time requirement manageable. 
All simulations reported here were performed on a single 
processor. 

The computation of the boundary integral in Eq. ([l} , 
requires the knowledge of the hydrodynamic traction 
jump across the interface (Af). This jump in traction 
is obtained from the membrane equilibrium condition, 
which requires that the hydrodynamic forces on any in- 
finitesimal area of the membrane be balanced by the elas- 
tic forces in the membrane [l^. Hence Af can be ob- 
tained from the knowledge of elastic stresses. The cap- 
sule membrane is modeled here as a neo-Hookean mem- 
brane with a shear modulus Gg [3. The elastic forces 
in the membrane are then computed using the approach 
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of Charrier et al. which is based on the principle 
of virtual work. Imp lementation details of this approach 
can be found in [ifl, [23| • 



C. Physical parameters 

We detail here the various physical parameters relevant 
for this study. To begin, we note that the rest shape of all 
capsules is spherical with radius a. This radius a will be 
taken as the unit of length throughout this work. Next, 
we introduce the undisturbed pressure driven flow, whose 
velocity in the flow direction x is given by: 



4{/o 
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(3) 



where C/q is the centerline velocity. Based on this flow 
profile, we obtain a wall shear rate of 7^, = AUq/H . Time 
in this work will be presented in units of 7~^. The shear 
modulus of the capsule will be expressed by the non- 
dimensional capillary number Ca = iijwO-/Gs, where fi 
is the viscosity of the suspending fluid. The capillary 
number can be viewed as the ratio of viscous and elas- 
tic stresses on the capsule. Since we study suspensions 
of binary mixtures of capsules with different rigidity, we 
remark that in a given flow, the values of Ca of the dif- 
ferent particles are not equal: particles with the lower Gs 
have a higher Ca and are termed floppy particles, while 
particles with the higher Gg have a lower Ca and are 
termed stiff particles. The ratio of the rigidities (Gs) of 
the stiff particle and the floppy particle will be denoted 
by R, with R > 1. Another important parameter charac- 
terizing a binary mixture is the number fraction of one of 
species in the suspension - in this work, we employ the 
number fraction of the floppy particles Xf. The number 
fraction of the stiff particles can be obtained from Xf as: 
Xg = 1 — Xf. The suspension as a whole is character- 
ized by its volume fraction 0, so for example Xf4) is the 
overall volume fraction of floppy particles. 

In the present work, we have neglected bending resis- 
tance of the capsule membranes. In such a case, com- 
pressive stresses in the membrane can lead to membrane 
buckling at both at low and high Ca [2l|. In the buck- 
led state, the numerical result is likely spurious and does 
not represent true physics. A simple way to alleviate this 
problem is to introduce an isotropic prestress in the mem- 
brane, which can also arise naturally in experiments due 
to osmotic effects [l^. As a result of the prestress, the 
rest radius of the capsule a is larger than the unstressed 
radius oq. The degree of inflation can be characterized 
by the inflation ratio beta = a/ao — 1. In pair collision 
studies, inflation ratios between 0.05—0.1 are typically re- 
quired to prevent buckling depending on the Ca, though, 
at higher Ca even these inflation ratios are insufficient 
[2]| . In this work, we set the inflation ratio to (3 — 0.1. 
We did not observe any buckling instability in our sim- 
ulations. The effect of preinflation on the cross-stream 
displacement in pair collisions is usually weak [2l|. We 



emphasize that the inflated radius a is taken as the rel- 
evant radius of the capsule; for example, Ca is defined 
based on the inflated radius a. 

Having defined all the important parameters, we next 
tabulate the set of simulations performed in this study 
with their specific choices of parameters (Table IJ). As 
can be seen, we have performed four main sets of simu- 
lations, denoted A, B, C and D in the table. These sets 
of simulations address the effect of the number fraction 
of the floppy particle Xf on the segregation behavior in 
combination with several other parameters including (j), 
Ca, R and the confinement ratio 2a/H. 



D. Simulation details 

We next outline the major steps involved in the solu- 
tion procedure. The first step involves the computation 
of the hydrodynamic traction jump Af - this is computed 
using the approach outlined above in Sec. (jll B[) . Follow- 
ing this, we compute the velocity at all surface element 
nodes employin g a discretized version of the boundary 
integral Eq. ([T]) . The surface element nodes are then 
advanced in time using the second order explicit Adams- 
Bashforth method. The time step At is set adaptively 
to satisfy: At = 0.5 Ca h, where h is the minimum node- 
to-node distance (the nearest nodes do not have to be 
on the same particle). Recall that time is expressed in 
units of 7^^, while distances are expressed in units of a. 
Although the above choice of time step At is sufficiently 
small in general, occasional particle-particle overlaps can 
occur. In principle the time step can be reduced further 
until overlaps are eliminated, but the computational cost 
can become exorbitant. Alternatively, the overlaps can 
be prevented by employing a short range repulsive force, 
though this is not very effective in three dimensional sim- 
ulations [11] . The best approach appears to be an overlap 
correction in an auxiliary step [235 . This approach is also 
frequently used in suspensions of rigid particles i25| . 
In the present work, like in suspensions of rigid particles 
[l^l, we not only correct the overlaps, but also maintain 
a minimum gap ((5mm) between the surfaces of two par- 
ticles. In the present case, we set the minimum gap to 
a small value of one percent of the (rest) particle diame- 
ter, i.e. Sniin — 0.02a. Maintaining the minimum gap is 
also justified in the absence of adaptive meshing, as that 
will be necessary to accurately resolve the flow dynamics 
in the small gaps. No adaptive meshing was employed 
in this work. The overlap correction procedure involves 
moving a pair of overlapping particles apart along their 
line of centers like a rigid particle until the minimum gap 
requirement is satisfied. Translating the capsules like a 
rigid particle in this auxiliary step has the benefit that 
the shapes of the particles remain unchanged, as is the 
orientation of the particles with respect to the flow. In 
general, multiple steps of the correction procedure are 
required, as the correction of overlap between one pair 
could result in other overlaps 2j|. This overlap correc- 
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TABLE I: Parameter specification in various simulation sets. 



Set 


Cai 


Ca.2 


R 


Xf 





2a/H 


Np 


A 


0.2 


0.5 


2.5 


- 1 


0.2 


0.197 


50 


B 


0.2 


0.5 


2.5 


- 1 


0.2 


0.3 


50 


C 


0.2 


0.5 


2.5 


- 1 


0.12 


0.197 


30 


D 


0.1 


0.6 


6.0 


- 1 


0.12 


0.197 


30 




(c) Xf = 0.8 



FIG. 2; (color online) Simulation snapshots in the shear plane 
{x — y) for different values of X/. The snapshots show ap- 
proximately 1.5 times the length of the simulation box in the 
X direction. The stiffer particles are light (blue online), while 
the floppy particles are dark (red online). 



tion step involves minimal displacement of the particles, 
on the order of 6min- Given that the maximum volume 
fraction studied in this work is (j) = 0.2, this procedure 
was rarely required - only a couple of minimum gap vi- 
olations were observed per time unit. Time integrations 
were performed for sufficiently long times (hundreds of 
time units) to allow development of an statistical steady 
state in the particle distribution. 



III. RESULTS 

We begin by discussing the results for simulation set 
A of table (|T|. In this set of runs, the stiffer particle has 
a capillary number of Ca = 0.2, while the floppy particle 



has a capillary number of Ca = 0.5, such that the rigid- 
ity ratio R = 2.5. The confinement ratio for this set of 
runs was fixed at 2a/ H = 0.197. The number fraction of 
floppy particles in the mixture {Xf) was varied between 
and 1 , with Xf = Q corresponding to a pure suspension 
of stiff particles (Ca = 0.2) and Xf — 1 corresponding to 
a pure suspension of floppy particles (Ca = 0.5). Several 
simulation snapshots for different values oi Xf are shown 
in Fig. ^ to help in visualizing the problem set up. We 
also remark that, throughout this article, we will often 
use the letter 'S' to denote stiff particles, and the letter 
'F' to denote floppy particles. 



A. Particle distribution in the wall normal 
direction 

The first suspension property we examine is the dis- 
tribution of particles in the wall normal direction y. To 
characterize this distribution of particles, we will employ 
three complementary measures. In the first measure, we 
compute the mean absolute distance of a species from 
the channel centerline, denoted by d; see Fig. ^^). This 
simple measure, though, does not provide detailed infor- 
mation about the actual distribution of particles along 
the wall normal direction. In order to quantify this dis- 
tribution, therefore, we will use two additional measures: 
the first is ^(y), which gives the volume fraction of a 
given species as a function of y-coordinate normalized 
by the mean volume fraction of that species, while the 
second measure is n{y), which gives the number density 
of the center of mass of a given species as a function of 
the y-coordinate normalized by the mean number den- 
sity n of that species. If the distribution were uniform, 
then both of these measures would have a value of unity 
throughout. Note that a subscript 's' or 'f may be used 
with or n for specifically referring to the stiff or the 
floppy particles, respectively; the same convention is fol- 
lowed later for denoting other quantities. We now briefly 
describe the procedure to compute n{y); a similar proce- 
dure is employed to compute (j){y). For this calculation, 
the channel height is flrst divided into bins and then the 
particles are assigned to bins based on their center of 
mass coordinates. The number of particles in each of the 
bins is then normalized by the value expected in that 
bin based on a uniform distribution in the wall normal 
direction. 
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(a) definition of d 
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1200 1600 2000 



(b) pure species {X f = 0, 1) 




400 800 1200 1600 2000 2400 

t 

(c) Xf = 0.2 

FIG. 3: (color online) (a) Definition of d: it gives the absolute 
distance from the centerline. (b) & (c) Time evolution of the 
mean absolute distance of a species from the centerline d non- 
dimensionalized by half the channel width H/2. In the plots 
'S' refers to the stiff particles (Ca = 0.2), while 'F' refers to 
the floppy particles (Ca — 0.5). In (b) results are shown for 
suspensions of pure species {Xf — &i Xf — 1), while in (c) 
results are shown for a suspension with Xf — 0.2. 



With the different measures for characterizing the par- 
ticle distribution specified, we turn to the results for these 
measures. We show the time evolution of d for the case 
of pure species (i.e. Xf = & Xf — 1) in Fig. ^jp) 
and for a mixture with Xf = 0.2 in Fig. It can 

be seen that in all cases d initially decreases with time, 
indicating that the particles are drifting toward the cen- 
terline. This behavior is not only consistent with results 
in suspensions of deformable particles subjected to pres- 
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FIG. 4: (color online) Steady state d as a function of Xf for 
the Ca = 0.2 and Ca — 0.5 particles. Also shown is the steady 
state d of the "overall" suspension. 



sure driven flows [l3|, l26|, but also with results in sus- 
pensions of rigid particles [13, HI]. However, as we shall 
see shortly, the behavior in suspensions of particle mix- 
tures is much richer. In suspensions of pure species (Fig. 
|31d), we find that d reaches an apparent steady state rel- 
atively quickly in about t — 500. In addition, the value 
of d for stiffer particles (Ca — 0.2) is found to be slightly 
higher than for floppy particles (Ca = 0.5), though the 
difference is not large. 

Next, we consider results for the mixture with Xf = 
0.2 (Fig. [S};). In this figure, wc show d for the stiff parti- 
cles (Ca = 0.2) and the floppy particles (Ca = 0.2) indi- 
vidually; in addition, d for the overall suspension is also 
shown. First of all, one can note that the d for Ca = 0.2 
particles and d for Ca = 0.5 particles are different at 
t = 0. This is due to the finite system size with the total 
number of particles being Np — 50, and the average being 
performed over 2 realizations; if large number of realiza- 
tions are used, then d of both the species will be the same, 
though this is not expected to affect the steady state dis- 
tribution. Returning to the figure, we observe that ini- 
tially both the stiff and the floppy particles drift toward 
the centerline. However, beyond approximately t = 600, 
in what appears to be a slowly evolving process, the stiff 
particles are displaced toward the walls by the floppy 
particles. The d of the overall suspension, though, shows 
much weaker evolution for t > 800, thereby confirming 
that this phase of the simulation mostly involves rear- 
rangement of positions between the two particle types. 
Munn and Dupin 13] saw a similar behavior in their sim- 
ulations on mixtures of deformable particles, which led 
them to the same conclusion as here. The mechanism 
underlying this observation is described in Sec. (|IV[) . 

We next summarize the results for the steady state d of 
both the species at all values of Xf studied here (Fig. [4]). 
Note that the steady state d was obtained by averaging 
over the last 1000 time (strain) units of the simulation. 
Two trends are immediately obvious in this plot. First is 
that d for both the stiff and the floppy particles increase 
with increasing Xf, though d for the whole suspension 
actually decreases with increasing Xf. The latter is ex- 
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pected as with increasing Xf, the suspension, on aver- 
age, becomes more deformable and hence has a lower d . 
An alternative interpretation of the above observation is 
that with increasing fraction of the floppy particles, the 
stiffer particles get displaced toward the wall, while with 
increasing fraction of the stiffer particles, the floppy par- 
ticles get displaced toward the centerline. The second 
observation is on the rate of increase or decrease in d of a 
given species with changing fraction of the other species: 
the rate of increase in d of stiffer particles is more rapid 
with increasing Xf than is the rate of decrease in d of 
floppy particles with decreasing Xf. This points toward 
a varying degree of segregation with changing Xf, which 
we will quantify shortly with an appropriate measure. 

To gain further insights into the particle distribution. 
Figs. ([U and ([5]) show fi{y) and (^(y), respectively. Note 
that in these figures, the channel wall is at the left and 
the centerline at the right - results in each half of the slit 
have been combined to enhance the statistics. We first 
examine Fig. ([5^), which shows n{y) for suspensions of 
pure species {Xf = Q and Xf = 1). In this plot two 
peaks are very apparent: a sharp peak close to the wall 
and a broad peak close to the centerline. Even though 
the peak close to the wall looks sharper in comparison to 
the one at the centerline, the width of this peak is rela- 
tively small; hence the volume fraction around the sharp 
peak is much lower in comparison to the volume fraction 
around the centerline (Fig. [5^). Also, the peak next to 
the wall is closer to it for the stiffer particles than for 
the floppy particles, thereby indicating that the cell free 
layer is wider in suspensions of more deformable particles 
consistent with other simulations [26|. The presence of 
peaks also points toward the formation of layers, though 
particle snapshots (see Fig. [5]) and animations reveal that 
only the peaks near the walls are representative of layers, 
while much of the interior region consists of disorganized 
particles. This is similar to observations is confined sus- 
pensions of drops [1^. The layering close to the walls in 
concentrated suspensions is easily explained by recalling 
that the wall is a no penetration boundary, thereby pro- 
moting the formation of layers next to it. It is also worth 
pointing out that the volume fraction at the centerline 
(Fig. is about twice the average volume fraction. 

Hence the suspension is concentrated around the center- 
line with the volume fraction there being (j) k, 0.4. 

We next examine the effect of varying Xf on h{y) and 
4>{y) of individual species. We focus first on the stiff par- 
ticles, for which fisiy) is shown in Fig. ^jp) and 4>s{y) in 
Fig. ([BJd). It is clearly seen that with increasing fraction 
of floppy particles, the concentration of stiff particles in- 
creases in the near wall region; at low values Xf this 
increase is mostly compensated by the decrease in con- 
centration in the region between the two peaks (one at 
the centerline and one next to the wall), while at higher 
values of X/, the stiffer particles are also slightly depleted 
in the region around the centerline. For the floppy par- 
ticles (Fig. [St and Ht), the exact opposite behavior is 
observed, namely, with increasing fraction of stiffer par- 
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FIG. 5: (color online) Number density distribution of a given 
species in the wall normal direction normalized by the mean 
number density of that species, n{y)- The channel wall is at 
y — and the centerline at j/ = 0.5H. Note that a uniform 
distribution will yield a value of one throughout, (a) Pure 
species: Xf = & Xf = 1, (b) ns{y) for stiH' particles at 
various values of Xf, and (c) hf{y) for floppy particles at 
various values of X/. 



tides (i.e. decreasing Xf), floppy particles get depleted 
in the region close to the walls and accumulate in the 
region between the two peaks; no significant changes in 
the concentration are observed in the region close to the 
centerline, except at very low values of X/. 

We turn now to quantifying the degree of segregation 
in the mixture. Here, we employ a particularly simple 
measure based on the volume fraction distribution of the 
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FIG. 6: (color online) Volume fraction distribution of a given 
species in the wall normal direction normalized by the mean 
volume fraction of that species, ^(y). The channel wall is at 
y = and the centerline at j/ = 0.5H. A uniform distribution 
will yield a value of one throughout, (a) Pure species: Xf — 
& Xf = 1, (b) 0s (j/) for stiff particles at various values of Xf, 
and (c) (j}f{y) for floppy particles at various values of Xf. 



two species as a function of the wall normal coordinate. 
To begin with, we compute Xf{y) defined as the frac- 
tion of the total particle volume contributed by floppy 
particles as a function of the wall normal coordinate: 



4)s{y) + (j)f{v) 



(4) 



This fraction normalized by the mean number fraction of 
the floppy particles Xf is denoted by Xf{y): Xf{y) = 
Xf{y)/Xf, which essentially gives the local enrichment 
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FIG. 7: (color online) (a) Local fraction of the total particle 
volume contributed by floppy particles normalized by its bulk 
value, Xf. (b) The enrichment factor E for both floppy (F) 
and stiff (S) particles as a function of Xf. 



(when Xf{y) > 1) or depletion (when Xf{y) < 1) of the 
floppy particles. We show in Fig. ([7^) the plot of Xf as 
a function of the wall normal coordinate (y) for various 
values of Xf. As expected, in all cases, floppy particles 
are enriched in the region close to the centerline, while 
the stiffer particles are enriched in the region close to the 
walls. If one's goal is to separate these two species, then 
an obvious way to achieve this is by separating the sus- 
pension present in these two regions. In this particular 
case, the plane separating these two regions is found to 
he y Ki H/4. We then define an enrichment factor E for 
each species as the ratio of its fractional contribution to 
the total particle volume in its enriched region normal- 
ized by the corresponding bulk value. For example the 
enrichment factor for floppy particles Ef is obtained as 



E 



f 



1 



f<j)f{y)dy \ 



(5) 



where the integrals in the above equation are restricted 
to the region where X/ > 1. We remark that, by defini- 
tion, E > 1. We show the plot of E for both species as a 
function of the mean fraction of the floppy particles Xf 
in Fig. ^h)- This result shows that the enrichment of 
either species increases as it becomes more dilute. More- 
over, the rate of this increase is much more rapid for the 
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FIG. 8: (color online) (a) Average velocity Vp of a particle of a 
given type in suspensions of pure species (S & F) as a function 
of y. Also shown are the velocities of an isolated (single) stiff 
(S*) and floppy (F*) particles in the same geometry and flow. 

(b) Relative change in the velocity of the stiff particles at a 
given Xf over the corresponding pure species result in (a). 

(c) Same as (b), but for the floppy particles. 



stiffer particles than for the floppy particles. This slow in- 
crease for the floppy particles is perhaps not too surpris- 
ing, as its enriched region is centered around the center- 
line, where even the stifFer particles have a propensity to 
accumulate, thereby reducing the relative enhancement 
of the floppy particles. 



B. Particle velocity in the wall normal direction 

So far we have focused on the particle distribution in 
the wall normal direction. It will also be of great inter- 
est to examine the particle velocity in the wall normal 
direction, as that is expected to play an important role 
in its distribution in the wall normal direction. We first 
show, in Fig. (|8K), the mean wall normal particle ve- 
locity Vp{y) in suspensions of pure species {Xj = and 
Xf = 1). This is obtained by first dividing the channel 
height into bins. A given particle's instantaneous wall 
normal velocity is then assigned to a bin based on the 
wall normal coordinate of its center of mass; the average 
velocity in a bin then gives the Vp in Fig. ([8^) . For com- 
parison, we also plot in the figure the quasi-steady wall 
normal velocity of an isolated capsule. This quasi-steady 
velocity is obtained by holding the center of mass of the 
isolated particle steady at a given location - at the end 
of each timestep, the particle was translated so that its 
center of mass returned to the original position. The re- 
sulting steady state particle velocity is referred to as the 
isolated particle quasi-steady migration velocity. 

The mean particle velocity Vp in the suspension ap- 
pears to show three different regimes (Fig. In the 
region close to the walls {y/H < 0.2), Vp is a decreas- 
ing function of the distance from the walls, much like 
the migration velocity of an isolated particle, though the 
velocity in the suspension is significantly lower. Recall 
that y/H — 0.2 approximately corresponds to the peak 
of h{y) and (f>{y) (Figs. [Hl&ISl)- In the intermediate re- 
gion {0.2H < y < 0.3-ff), Vp is approximately a constant, 
and can even be a slightly increasing function of dis- 
tance from the wall?/. In the region close to the centerline 
{0.3H < y < 0.5H), the velocity decreases linearly with 
distance. This behavior is similar to that of an isolated 
particle, whose migration velocity close to the centerline 
also decreases linearly with distance. This results from 
the fact the shear rate decreases linearly with distance in 
Poiseuille flow, while the gradient of the shear rate is a 
constant. The migration velocity, which, to leading or- 
der near the centerline, is the product of these two terms, 
therefore decreases linearly with distance [s^l- Perhaps 
the most interesting feature in Fig. (|8^) is the fact the 
mean particle velocity in the suspension can be higher 
than the corresponding isolated particle migration veloc- 
ity. This effect certainly arises due to particle-particle 
interactions and is discussed further below. 

We next turn to the effect of AT/ on Vp of particles of 
a given species in the suspension. To quantify this, we 
compute the relative change in the velocity Avp/vr = 
(vp — Vr)/vr, where u,- refers to Vp in the corresponding 
pure species suspension. We show plots of Avp/vr for the 
stiff particles at several values of Xf in Fig. ^jp), while 
the same plots for the floppy particles are shown in Fig. 
dlt). Focusing first on the stiff particles, we find that 
their mean wall normal velocity Vp is greatly reduced, in 
some cases by nearly 50%, as the fraction of floppy par- 
ticles in the suspension increases. The most significant 
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FIG. 9: (color online) Probability distribution of the wall 
normal particle velocity Vp at several values of y/H . (a) pure 
stiff particles {Xf = 0), and (b) pure floppy particles (Xf — 
1). 



decrease in the velocity is observed in the second region 
discussed above, where Vp of the particles in suspensions 
of pure species is approximately a constant (Fig. 
A significant decrease is also observed in a thin region 
around the centerline, though, Vp itself is very small here. 
To gain some insights into this behavior, we remark that 
with increasing fraction of the floppy particles, a stiff par- 
ticle is expected to "collide" more frequently with floppy 
particles than with particles of the same kind, particu- 
larly from its side facing the centerline where most of the 
floppy particles are concentrated. Therefore, it appears 
that the mean additional effect of these heterogeneous 
collisions is to push the stiffer particles toward the wall. 
Later, in Sec. (|IVp . it will be shown that this observation 
is consistent with simple model studies on homogeneous 
and heterogeneous pair collisions. We next consider the 
effect of Xf on Vp of floppy particles (Fig. [8j;). It is 
immediately obvious that, in comparison to Vp of stiffer 
particles, Vp of floppy particles is a weak function of Xf 
and no particular conclusion can be drawn based on these 
trends. 

Above, we examined the mean particle velocity Vp. 
This particle velocity also exhibits considerable fluctu- 
ations about this mean. To emphasize this, we show in 
Fig. ^ the probability distribution of particle velocities 
P{vp) on a log-linear scale at various wall normal coor- 
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dinates in suspensions of pure stiff and floppy particles, 
i.e. Xf = and Xf — 1, respectively. Not only do these 
plots confirm that there are considerable fluctuations in 
particle velocity (this is quantified below), but the exact 
nature of these distributions is interesting in itself and 
deserves a brief discussion. Perhaps the most interesting 
observation in these plots is that the distribution of posi- 
tive velocities for particles close to the wall {y = 0.2H in 
the figure) is non-Gaussian, and is, in fact, much closer 
to an exponential distribution. At larger values of y/H, 
this changes to a Gaussian distribution. In contrast, the 
distribution of negative velocities appears to be Gaussian 
at all values of y/H. Since one usually expects a Gaus- 
sian distribution for a variable that is obtained as a sum 
of many uncorrelated parts (in view of the central limit 
theorem [31]), the appearance of an exponential distribu- 
tion perhaps implies contributions from correlated parts 
[3^. In the present case, the underlying physics of the 
problem provides a basis for these observations. We note 
that a particle is expected to show a non-zero wall normal 
velocity both due to "random" collisions (interactions) 
with its neighbors, as well as due to the non-random ef- 
fects of the shear rate gradient and hydrodynamic wall 
repulsion [s^, HI] . In light of the above discussion, there- 
fore, the former contribution can be expected to lead to a 
Gaussian distribution, while the latter to a non-Gaussian 
one. Since both the shear rate gradient and wall repul- 
sion lead to positive velocities, it is therefore not surpris- 
ing that one observes a non-Gaussian positive velocity 
distribution, especially close to the walls where the non- 
random effect is dominant. Negative velocities, on the 
other hand, can result only due to interactions of a par- 
ticle with its neighbors, and therefore can be expected to 
have a Gaussian distribution. 

We next quantify the fluctuation in the particle veloc- 
ity by computing its standard deviation, denoted by v'p. 
We flrst show in Fig. pI7b ) plots of v'p as a function of y in 
suspensions of pure species. As can be seen, the fluctua- 
tions in the velocity are much larger than the correspond- 
ing mean, especially in the region close to the centerline. 
This is particularly true for the stiff particles as they not 
only show a smaller Vp, but also a higher v'^. We also 
find in the plot that the fluctuations show a minimum in 
the region near the wall (y « 0.2H). This minimum in 
fluctuation appears to result from the orderly motion of 
particles in the layer formed next to the wall (Fig. |6^). 
Moving further away from the wall, we observe a maxi- 
mum in the fluctuation v'^ at approximately y = O.SH. 

This position corresponds well with the minimum in 4'{y) 
in Fig. (|6]). The maximum in the fluctuation implies 
that a particle at this position undergoes more frequent 
collisions with other particles than do particles at neigh- 
boring positions, and hence a particle in this region is 
expected to drift away from this region. It is not sur- 
prising, therefore, that the position of the maximum in 
velocity fluctuation corresponds well with the minimum 
in 0(y). In addition, we observe that the fluctuation am- 
plitude in the region close to the centerline is non-zero. 
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FIG. 10: (color online) Fluctuation (standard deviation) of 
the particle velocity of a given species as a function oi y/H. 
(a) In suspensions of pure species {Xf = & Xf = 1), (b) for 
stiffer particles at various values of Xf, (c) for floppy particles 
at various values of X/. 



stiffer particles was always higher than the corresponding 
fluctuations in the velocity of floppy particles. 

Based on the particle velocity fluctuations above, 
which in turn are rooted in particle collisions, one may 
explain the enhancement in Vp of a particle in the sus- 
pension over the corresponding isolated particle near the 
centerline (Fig. Because of the variation in the ve- 
locity fluctuations (collisions) as a function of y, one can 
expect a particle to drift down the gradient of velocity 
fluctuations (collision) [l^, 113]; the velocity fluctuation 
profile in Fig. ([TOk ) then implies a drift toward the cen- 
terline due to this effect, thereby leading to an enhance- 
ment in the Vp. 

Results with other simulation parameters in Tb. 
show a qualitatively similar behavior and are presented in 
the appendix for completeness. In the following section, 
we explore the mechanism responsible for the segregation 
behavior. 



IV. DEVELOPMENT OF MECHANISTIC 
MODEL 

In this section, we elucidate the mechanism that leads 
to the segregation between the stiff and the floppy parti- 
cles in suspensions of mixtures. In particular, we seek to 
explain the following two observations on the segregation 
behavior: (i) in a suspension of primarily floppy particles 
(i.e. large Xf), the stiff particles accumulate in the near 
wall region, while being depleted in the region around 
the centerline, and (ii) in a suspension of primarily stiff 
particles, the floppy particles are depleted from the near 
wall region accompanied by an increased accumulation in 
the region near the centerline. In these statements, the 
degree of depletion and excess of a given species in a par- 
ticular region is with respect to the distribution of that 
species in its pure suspension (i.e. Xf = or Xf — 1). 
We will show that these observations can be qualitatively 
explained with a unified mechanism incorporating the 
wall-induced migration of deformable suspended parti- 
cles away from the wall and the particle fluxes associated 
with heterogeneous and homogeneous pair collisions. 



Framework 



and stays flnite even at the centerline {y = 0.5H) where 
the local shear rate is zero. This is a well known result 
in monodisperse suspensions and has been attributed to 
non-local effects due to the finite particle size [ll|, [2^ . 

We next summarize the effect oi Xf on the velocity 
fiuctuations of the stiffer particles in Fig. ([TUb). and for 
the floppy particles in Fig. ([TUb ). These plots show that 
the velocity fluctuation is a weak function oi Xf, though, 
in general, the fluctuations were higher in mixtures than 
in the corresponding suspensions of pure species. We 
also point out that the fluctuations in the velocity of 



To make progress toward elucidating the segregation 
mechanism, it will be helpful to flrst develop a model- 
ing framework with which to analyze and interpret the 
simulations results. In particular, we focus on the parti- 
cle motion in the wall normal direction and identify the 
distinct mechanisms that result in such a motion. 

An isolated deformable particle in a pressure driven 
flow acquires a non-zero wall normal migration velocity 
both due to hydrodynamic interactions with the walls 
(wall repulsion) 33] and due to the curvature in the flow 
proflle (shear rate gradient) [s^]- In Sec. (jIII[) . we pre- 
sented the quasi-steady migration velocity of an isolated 
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capsule in Poiseuille flow; see Fig. (jS]). For the present 
section, we denote this velocity by , where the su- 
perscript 'S' signifies the origin as the single (isolated) 
particle motion. 

In a suspension of particles, the interparticle interac- 
tions also contribute to the particle's wall normal drift 
velocity (29j . In general, developing exact expressions 
for the particle velocity in the suspension is difficult, 
although for dilute solutions the treatment is relatively 
simple. This will prove sufficient for our present, qualita- 
tive purposes. In the dilute limit, as an approximation, 
we can sum velocity contributions resulting from the iso- 
lated particle motion and due to contributions from par- 
ticle interactions treated in a pairwise fashion. The ne- 
glect of three-particle and higher order interactions will 
cause an error of 0(0^), which is small in dilute suspen- 
sions. We will denote the contribution to a particle's 
velocity resulting from pair collisions by , where the 
superscript 'P' signifies the origin of the term being pair 
collisions. We can combine the contributions from pair 
collisions and single particle motion to obtain the total 
wall normal particle drift velocity as: 



V ^ V" 



(6) 



We next develop an expression for the velocity of a par- 
ticle resulting due to pair collisions (V^). 

The (ensemble-averaged) velocity of a particle aris- 
ing due to pair collisions in a dilute suspension can be 
expressed by a precise integral expression; see, e.g., Li 
and Pozrikidis |29|. For our qualitative discussion in the 
present work, however, we employ and expand on the 
much simpler phenomenological model of Phillips et al. 
27| . In this model, the number of collisions experienced 
by a test particle is taken to scale as 71/). The variation 
of this collision frequency over a distance of 0{a) is then 
given by a\7{'y(f)), where V here reduces to the derivative 
with respect to y. Furthermore, each of these collisions 
are assumed to cause a displacement of 0(a), thereby 
leading to a particle velocity given by: 



(7) 



where X is a dimensionless constant that may depend on 
the properties of the particles in the suspension (e.g. cap- 
illary number, shape). A few words are in order about 
the interpretation of this parameter. We note that the 
above expression was derived under the assumption that 
the particle undergoes a displacement of 0{a) in pair 
collisions, which is an appropriate order of magnitude 
estimate, though it does not provide precise informa- 
tion about the actual characteristic displacement in pair 
collisions. This information takes additional significance 
when considering a suspension of mixtures, as the char- 
acteristic displacement in pair collisions is expected to 
be different when a test particle undergoes collision with 
a particle of the same type (homogeneous collision) or 
with a particle of the other type (heterogeneous collision). 
Furthermore, the characteristic displacement will also de- 
pend on the identity of the test particle, i.e. whether it 



is stiff or floppy. In order to make the above model (Eq. 
[7]) cognizant of this information, we view 'Ka' as the 
characteristic displacement of a particle in pair collisions. 
For example, when considering pair collisions between 
smooth rigid spheres, we have K — a,s the characteris- 
tic displacement in such collisions vanishes owing to the 
symmetry of the geometry and the Stokes flow reversibil- 
ity. With this aspect clarified, we generalize Eq. ([7]) for 
a suspension of mixtures to obtain: 



Vf = -a' (if..V(7<^s) + if./V (j^f)) 



(8a) 



^a^{KfsVij(l,s) + KffVij(l,f)), (8b) 



where the subscript 's' and 'f denote the stiff and the 
floppy particle, respectively. We next discuss results from 
several pair collision studies as they provide important 
insights into the relative values of the characteristic dis- 
placements in various scenarios, i.e. it provides informa- 
tion about the relative values of the various constants K 
in equation ([8]). 



B. Pair collisions 

We presented above the particle velocity resulting from 
pair collisions between particles, which obviously plays 
an important role in the particle distribution along the 
wall normal direction. To gain more insights into this 
process, it will be helpful to examine some pair collision 
results, particularly for pair collisions between two differ- 
ent species, as no prior study appears to have addressed 
this problem. We begin by describing the problem setup 
for pair collision studies. In all studies, a large cubic 
simulation box of side 30a was considered, and the par- 
ticles were placed on either side of the central plane with 
a given initial separation. Further, in all these studies, 
only simple shear fiows were considered. These choices 
for the problem setup minimize the single particle migra- 
tion effect, and hence provide insights into the effect of 
particle-particle interactions (collisions) . 

We begin by considering pair collisions between the 
same species (homogeneous collisions) and show in Fig. 
PTk ) the displacement of the 'top' particle in the wall 
normal direction relative to its initial position (y — yi) 
as a function of time. For this study and those below, 
we set the initial separation between the two particles in 
the wall normal direction as Ayi = 0.5a. Other studies 
with Ayi — 0.25a and Ayi = 0.75a gave similar trends. 
In the above figure, it is clearly seen that in this range of 
Ca (0.1 < Ca < 0.6), the final displacement {Sea = yoo — 
yi) increases with decreasing Ca. We remark that 5 00 is 
known to be a non-monotonic function of Ca, such that 
at low enough enough Ca, 5oo decreases with decreasing 
Ca |20i]. 

We next consider collisions between two particles with 
different Ca (Fig. [Tib). From these results we see that 
the stiffer particle undergoes much larger displacement 
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FIG. 11: (color online) Pair collisions: Displacement of a par- 
ticle relative to its initial position (y — j/i) with time in (a) 
homogeneous pair collisions and (b) heterogeneous pair colli- 
sions. Wall normal velocity of a particle as a function of time 
in (c) homogeneous pair collisions and (d) heterogeneous pair 
collisions. 



(Joo) than the corresponding floppy particle involved in 
the collision. In addition, the difference in the final dis- 



placements of the two species is larger for larger rigidity 
ratio system. Further insights are obtained by comparing 
the instantaneous wall normal velocities of the particles 
during the collision event for the cases in Figs. pTk ) 
and pTb ). which are shown in Figs. ([Tib) and ([TTJi), 
respectively. Focusing first on the result for particle ve- 
locity in homogeneous pair collisions (Fig. ITTb). we find 
the particle velocity is positive in the approach part of 
the collision trajectory, while it is negative in the reced- 
ing part of the collision trajectory. As a result of this 
feature, a fraction of the maximum displacement in Fig. 
([TTb) is recovered in the receding part of the collision tra- 
jectory. It is obvious from Fig. (fTTb ') that the magnitude 
of velocity in either part of the trajectory is higher for 
the stiffer particle than the fioppy particle in homoge- 
neous collisions. In contrast, in heterogeneous collisions 
(Fig. [TTJi), while the magnitude of the velocity is higher 
for the stiff particle in the approach part of the trajec- 
tory, in the receding part it is the floppy particle which 
exhibits a greater magnitude of velocity. This reversal in 
the trends of the velocity in the receding part of the tra- 
jectory leads to the observed behavior in heterogeneous 
pair collisions. Summarizing the results of this section, 
we have the following relationships between the various 
constants K in Eq. ([S|): (i) K^f > Kfg, (n) K^f > K^s, 
and (iii) Kfs < Kff. We next discuss the implications of 
these pair collision results on the segregation behavior. 



C. Mechanism of flow induced segregation 

Based on the results of pair collision studies and the 
single particle migration velocity, we are now in a position 
to identify the mechanism responsible for the segregation 
behavior. In a confined suspension, the wall plays an im- 
portant role as it breaks the translational invariance in 
the wall normal direction. So, we divide the overall sys- 
tem domain in two parts: first is the near wall region, 
and second is the interior region. These two regions are 
pictorially shown in Fig. (|12p . where the near wall region 
is labeled as 'A', while the interior region is labeled as 'B'. 
We first focus on the near wall region as that turns out 
to be critical to the segregation behavior. For a particle 
in the near wall region A (e.g., shaded particles in Fig. 
[T^ , the net of effect of pair collisions is to push to parti- 
cle toward the wall, as — V^jcf)) is negative in this region. 
In other words, a particle in this region undergoes colli- 
sions with other particles present only on its upper half 
and consequently gets pushed toward the wall. There- 
fore, the net effect of these pair collisions is to prevent 
the 'escape' of the particles in the near wall region A 
to the interior region B. In contrast, the single particle 
velocity , which is always toward the centerline, aids 
the particle in 'escaping' to the interior region B. Since 
the single particle migration velocity depends only 
on the characteristics of the particle, it is independent 
oi Xf. On the other hand, , the velocity arising due 
to pair collisions, depends strongly on Xf in view of the 
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FIG. 12: (color online) A cartoon depicting the characteristics 
of the neighborhood of a particle of given type (shaded) in the 
near wall region at various values of X/. The stiff particle is 
shown light (blue online), while the floppy particle is shown 
dark (red online). 



results of homogeneous and heterogeneous pair collisions 
discussed above, i.e. due to differences in values of Kss, 
Kff, Kgf and Kfg. Therefore, the frequency with which 
a particle present in the region A will escape to the region 
B will depend on X/ primarily due to the dependence of 
on Xf. 

To aid in visualizing the problem at various values of 
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FIG. 13: (color orfline) (a) A schematic defining the near wall 
region A and the interior region B. (b) Transition probability 
/ab defined as the fraction of particles starting in region A 
crossing over to region Hint = 100. /ab is shown as function 
of Xf for both the stiff and the fioppy particles. 



Xf, we consider four cases as shown in Fig. ((T2|) : (a) sus- 
pension of stiff particles only (Xf — 0), (b) suspension 
that is dilute in stiff particles (Xf — > 1), (c) suspension 
of floppy particles only {X f = 1), (d) suspension that 
is dilute in floppy particles {Xf — > 0). In the first two 
cases, our focus is on the behavior of the stiff particle, 
while in the last two cases, our focus is on the behavior 
of the floppy particle; in all cases the particle of interest 
is shaded. Focusing first on a stiff particle, its velocity 

in cases (a) and (b) is given by = —Kssa'^'^{i(f)) 
and = -Ksfa'^Vijcj)), respectively. In Sec. (IIV^ . 
we concluded that Kss < Ksf, and hence that the stiff 
particle in case (a) faces a smaller barrier in escaping to 
region B than the stiff particle in case (b). Therefore, we 
conclude that a stiff particle gets increasingly localized 
in the near wall region as the fraction of floppy particles 
in the suspension increases. We next turn our attention 
to the behavior of the floppy particle and consider cases 
(c) and (d) in Fig. p^ . The velocity of the floppy par- 
ticle in these cases is given by — —KffO?V{'^<l)) and 

= —Kfsa?'^{"f(j)), respectively. The results of Sec. 
(lIVBp indicate that Kff > Kfs, and consequently that 
the frequency of escape of the floppy particle to region 
B from region A should increase as the fraction of stiff 
particles in the suspension increases. 

To verify the mechanism presented above, we compute 
the probability of transition of a particle of a given type 
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from region A to region B as a function oi Xf (Fig. [TSbl. 
More precisely, this transition probability, denoted by 
Jab , gives the fraction of particles that are present in re- 
gion A at a given time, but cross over to region B within 
100 time units from that time. We have defined the near 
wall region A as the region extending to a distance of 
0.3H from either walls, while the remainder of the do- 
main is defined as region B (Fig. [T5kV The choice of 
0.3H comes from the minimum of (j) in Fig. ([6^) and h 
in Fig. ([5^), which essentially separates the near wall re- 
gion and the interior region. The transition probability 
/ab in Fig. ([T3b ) clearly clearly shows trends consis- 
tent with the mechanism described above. We focus first 
on the trends for stiff particles. At Xf — 0, we find 
/1b ~ 0.081, which gradually decreases to /^^ « 0.015 
at Xf = 0.8, which amounts to a factor of 5.4 reduction 
in the transition probability. This result clearly suggests 
that stiff particles get localized in the near wall region as 
the fraction of floppy particles increases, which must be 
due to the effect of pair collisions, as the isolated particle 
migration velocity is independent oi Xf. As discussed 
earlier, this is the expected trend based on the relative 
values of Kgs and Ksf, as Ksf > Kss- We next focus on 
the trends of Jab for the floppy particle in Fig. (fT3b). At 
Xf = 1, we find Jab ~ 0.083, which gradually increases 
to « 0.1413 at Xf = 0.2, which amounts to a factor 
of 1.7 increase in the transition probability. This implies 
that, with increasing fraction of the stiff particles, the 
floppy particles escape to the interior region more eas- 
ily. This again is expected based on the results of pair 
collisions as Kfs < Kff. 

It is also interesting to note that the Jab for both the 
stiff and the floppy particles are very close in their cor- 
responding suspensions of pure species, i.e. Xf — Q and 
Xf = \, respectively. Results of pair collisions above in- 
dicate that the Kss and Kff are about the same, though 
the single particle migration velocity is larger for the 
floppy particle. A larger for the floppy particle, in 
principle, should lead to a higher Jab for that particle. 
On the other hand, a larger may lead to more frequent 
collisions with a relatively smaller Ay^, the initial sepa- 
ration between the two colliding particles in the wall nor- 
mal direction, as this velocity will cause two particles to 
approach each other in the wall normal direction. Since 
collisions with smaller /S.V j lead to larger cross stream 
displacement Saci{Ayi) [21|, the effect of larger may 
be nullifled by the effect of larger Sooi^Ui)- This effect 
enters our model through the effect of Sooi^Vi) on the 
constant K, as the constant K is essentially the weighted 
average of Soo{^yi), where the weight is the probability 
for a collision with a given Ayi to occur. Computing this 
weighted average rigorously is beyond the scope of the 
present work. 

In addition to the mechanism outlined above, we must 
also address the effects of excluded volume, which is only 
expected to amplify the effects of this mechanism. This 
is because if the floppy particles find it easier to escape to 
the interior than the stiff particles, then their enhanced 



presence in the interior will only make the escape of the 
stiff particles to the interior more difficult, as only a lim- 
ited volume fraction of particles can be present in the 
interior. So far we have been treating the motion of stiff 
particles and floppy particles separately. The excluded 
volume effect essentially couples the motion of these two 
species, which is intuitively expected to be the case in a 
suspension of particles, especially at non-dilute volume 
fractions. 

So far we have focused entirely on the near wall region. 
We now briefly discuss the dynamics in the interior re- 
gion. In this region, a particle has near neighbors both 
in its upper and lower half, and hence the gradient in the 
collision frequency aV(70) is expected to be small rela- 
tive to the collision frequency 70 itself (the latter collision 
frequency is only an approximation as it does not account 
for flnite particle size; see Sec. IIII B[) . When coupled with 
the fact that the single particle velocity is also small far 
away from the walls (Sec. IIII B|) . we conclude that the 
particle motion in this region is expected to be dominated 
by shear-induced diffusion. Therefore, the characteristics 
of the particle motion in this region are not expected to 
play a dominant role in the segregation between the two 
species. Nonetheless, the diffusive behavior in this region 
is necessary for transporting the particles from the bulk 
to the near wall region, and hence indirectly plays a role 
in the observed segregation between the two species. 

We conclude this section by addressing an interesting 
feature concerning the time evolution of the mean abso- 
lute distance d of a species in a mixture. In Sec. (|III Ap . 
we observed that in a mixture of stiff and floppy particles, 
d of both species initially decreases with time (Fig. [3];). 
But, beyond approximately t = 600 in Fig. the Jof 
stiff particles was found to gradually increase with time. 
This phase of the process was termed as a reorganiza- 
tion phase between the stiff and the floppy particles, as 
d of the overall suspension was essentially steady in this 
phase of the simulation. Based on the mechanism out- 
lined above, this reorganization phase can be associated 
with the transport of the stiff particles in the interior re- 
gion to the near wall region, where they are expected to 
become nearly irreversibly trapped. Therefore, the time 
scale after which this reorganization phase is noticeable 
in simulations can be associated with the time required 
for a particle near the centerline to diffuse to the near 
wall region. We therefore see that the mechanism out- 
lined here is capable of addressing even the finer aspects 
of the segregation behavior. 



V. CONCLUSIONS 

In this article we investigated the segregation behav- 
ior in suspensions of binary mixtures of Neo-Hookean 
capsules subjected to pressure driven fiow in a planar 
slit. The two species in the binary mixture have different 
membrane rigidities as characterized by their separate 
capillary numbers. Detailed boundary integral Simula- 
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tions revealed that with increasing fraction of the floppy 
particles in the suspension (i.e. with increasing Xf), there 
is an enhancement in the concentration of the stiff par- 
ticles in the near wall region. On the other hand, as the 
fraction of stiff particles in the suspension increases (i.e. 
Xf decreases), a depletion in the concentration of floppy 
particles in the near wall region is observed. As a result 
of this behavior, stiff particles segregate in the near wall 
region in a suspension of primarily floppy particles, while 
floppy particles segregate around the centerline in a sus- 
pension of primarily stiff particles. A novel and detailed 
mechanism based on the degree of localization of a par- 
ticle of a given type in the near wall region as a function 
of Xf was developed. The three main ingredients of this 
mechanism which governs the localization of a particle 
in the near wall region are: (i) the single particle migra- 
tion velocity toward the centerline, which provides an av- 
enue for the particle to escape to the interior, (ii) the re- 
peated pair collisions of a particle in the near wall region 
with its neighboring particles, which pushes it toward 
the wall, thereby preventing its escape to the interior 
and (iii) the difference in final displacement between ho- 
mogeneous and heterogeneous pair collisions. The single 
particle migration velocity is independent of X/, while 
the cross-stream displacement depends strongly on Xf. 
In particular, results of heterogenous pair collisions in- 
volving a stiff and a floppy particle showed that the stiff 
particle experiences a substantially larger cross-stream 
displacement as a result of the collision than the floppy 
particle. This result indicates, consistent with simulation 
results, that as the fraction of floppy particles in the sus- 
pension increases, a stiff particle in the near wall region 
will get increasingly localized there due to the strong ef- 
fect of heterogeneous pair collisions on its wall-normal 
displacement. On the other hand, as the fraction of stiff 
particles in the suspension increases, a floppy particle in 
the near wall region can escape to the interior much more 
easily due to the much weaker effect of heterogenous pair 
collisions on its wall-normal displacement. This mecha- 
nism, therefore, consistently explains the features of the 
segregation behavior observed in detailed numerical sim- 
ulations. Important topics for future work will be the 
extent to which the same picture can explain segregation 
by size and shape as well as segregation in the (nomi- 
nal) absence of shear rate gradients, as in plane Couette 
flow. It will also be of interest to formalize the quali- 
tative picture into a stochastic process model that more 
quantitatively captures the simulation and experimental 
results. 
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FIG. 14: (color online) (a) Steady state d as a function of 
Xf for stiff and floppy particles at two different confinement 
ratios. Enrichment factor as a function of Xf for (b) stiff 
particles, and (c) floppy particles. 



Appendix: Effect of confinement ratio, volume 
fraction and rigidity ratio 

We present here some additional results of simulations 
with different confinement ratio, volume fraction or rigid- 
ity ratio than were used in the main text. The results are 
all consistent with those results and are included here for 
completeness. 

To examine the effect of confinement ratio on the 
results, we consider suspensions with a volume frac- 
tion oi (f) — 0.2 and two different confinement ratios: 
2a/H = 0.197 and 2a/H = 0.3. The rest of the pa- 
rameters for these two simulation sets are listed in Tb 
(HI as sets A and B. We first show the effect of confine- 
ment ratio on d for both the stiff and the floppy particles 
in Fig. as a function oi Xf. As is obvious, d of 

both the species decreases with increasing confinement 
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FIG. 15: (color online) (a) Steady state d as a function of Xf 
for stiff and floppy particles at two different volume fractions. 
Enrichment factor as a function of Xf for (b) stiff particles, 
and (c) floppy particles. 
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FIG. 16: (color online) (a) Steady state d as a function of 
Xf for stiff and floppy particles as noted in the figure key. 
Enrichment factor as a function of X/ for (b) stiff particles, 
and (c) floppy particles. 



ratio, i.e. an increasing fraction of the particles accumu- 
late near the centerline with increasing 2a/ H. Moreover, 
the trend in d of either species with varying Xf is similar 
in both cases - with increasing fraction of floppy parti- 
cles (increasing Xf), the stiffcr particle gets increasingly 
displaced toward the wall, while with increasing fraction 
of stifFer particles (decreasing Xf), the floppy particles 
get increasingly displaced toward the centerline. 

For a more quantitative measure of the segregation be- 
havior, we compute the enrichment Es and Ef of the 
two species as a function of Xf (Figs. [Tib. fT4b). These 
plots show that E^ increases more rapidly with Xf in 
more confined suspension {2a/ H = 0.3), while Ef in- 
creases more rapidly with decreasing ^/ in less confined 
suspension {2a /H — 0.197). In other words, based on 
this measure, more confined suspensions lead to a higher 
maximum segregation of the stiff particles, while less con- 
fined suspensions leads to higher maximum segregation 



of floppy particles. 

We now we briefly examine the effect of volume fraction 
on the segregation behavior in the suspension. Results 
will be reported here for a suspension with a volume frac- 
tion of = 0.12 and compared with the results oi 4> = 0.2 
suspension discussed above in Sec. (IIIII) . The parameters 
for this simulation set are listed in Tb. (jH in the row 
labeled "C", which are the same as in set A except for 
the volume fraction. 

We begin by showing d in Fig. (fTSk ) for both the stiff 
(Ca — 0.2) and the floppy (Ca = 0.5) particles in the 
suspension as a function of Xf. For comparison, we also 
show d for the (f> = 0.2 suspension. At any given vol- 
ume fraction, d for stiffer particles is found to be higher 
than for floppy particles. Moreover, as the fraction of 
the floppy particles Xf increases, d for stiffer particles 
increases, i.e. they get increasingly displaced toward the 
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wall. On the other hand, as the fraction of stiff particles 
in the suspension increases {Xf decreases), J for floppy 
particles decreases, i.e. they get displaced toward the 
centerline. All these trends in cj) = 0.12 suspensions are 
similar to the trends in (p = 0.2 suspensions. It is also 
obvious in the plot that at any given value of Xf, the 
d for a species increases with increasing 0. This is an ex- 
pected result, as at higher volume fractions the increased 
particle-particle interactions reduces the fraction of par- 
ticles accumulating near the centerline, thereby leading 
to an increase in d . 

To quantify the degree of segregation, we compute the 
enrichment E for both the species as a function oi Xf 
(Figs. ITSbfTSb). These plots confirm that the maxi- 
mum enrichment of a given species is observed when that 
species is present in dilute amounts. Furthermore, the 
maximum enrichment of the stiff particles is higher for 
the (j) — 0.2 suspension, while the maximum enrichment 
of floppy particles is higher for the </> = 0.12 suspension. 

Finally, we consider the effect of Ca as well as that of 
the rigidity ratio R on the segregation behavior. For this 
section, we compare the results of set C and set D in Tb. 
(U]) - the volume fraction in both the simulations were 



(p = 0.12, while the stiff particles have a Ca of 0.2 and 
0.1, and the floppy particles have a Ca of 0.5 and 0.6, 
respectively. These sets of Ca for the stiff and the floppy 
particles give a rigidity ratio of i? = 2.5 and i? = 6, 
respectively. 

We begin by showing d for both the stiff and floppy 
particles in these simulations in Fig. ([TB^). One can 
observe in the plot that d for floppy particles is approx- 
imately the same in these two simulations. In contrast, 
d for the stiffer particles in the two cases show a marked 
difference as a function oi Xf - at higher values of R, 
the stiffer particle gets increasingly displaced toward the 
walls as Xf increases. For quantifying the degree of seg- 
regation between the two species, we again compute the 
enrichment factor E for both the species as function of 
Xf (Figs. [TBb . [TBb). As can be observed. Eg shows a 
much larger increase with increasing Xf for R — 6 case 
in comparison to the R = 2.5 case. In contrast, the en- 
richment of the floppy particles, Ef, is found to be nearly 
the same in both cases. These results suggest that the 
degree of segregation of the stiffer particle increases as 
the rigidity ratio increases, particularly at higher values 
oiXf. 
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